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1.  INTRODUCTION 


Tliere  are  now  three  basic  approaches  to  finite  element  analvsis.  In 
ail  approaches  the  domain  is  divided  into  simple  convex  subdomains  (usuallv 
triangles  or  rectangles  in  two  dimensions,  and  tetrahedra  or  bricks  in  three 
dimensions)  and  over  each  subdomain  Che  unknown  is  approximated  bv  a  (local) 
basis  function  (usually  a  polynomial  of  degree  p) .  Basis  functions  are 
required  to  meet  continuouslv  at  boundaries  of  subdomains  in  the  case  of  planar 
or  3  dimensional  elasticity,  or  smoothly  in  the  case  of  plate  bending.  The 
approaches  are: 

1.  The  h-version  of  the  finite  element  method.  In  this  approach  the  degree  p 
of  Che  approximating  polynomial  is  kept  fixed,  usually  at  some  low  number 
such  as  2  or  3.  Convergence  is  achieved  by  allowing  h,  the  maximum  diameter 
of  the  convex  subdomains,  to  go  to  zero.  Estimates  for  the  error  in  energy 
have  long  been  known  [1,  2,  3].  In  all  of  these  estimates  p  is  assumed  to 
be  fixed  and  the  error  estimate  is  as%-mptotic  in  h,  as  h  goes  to  zero. 

2.  The  p-version  of  the  finite  element  method.  In  this  approach  the  subdivi¬ 

sion  of  the  domain  Cl  is  kept  fixed  but  p  is  allowed  to  increase  until  a 
desired  accuracy  is  attained.  The  p-version  is  reminiscent  of  the  Ritz 
method  for  solving  partial  differential  equations  but  with  a  crucial  distinc¬ 
tion  between  the  two  methods.  In  the  Ritz  method  a  single  polynomial 
approximation  is  used  over  the  entire  domain  12  (12,  in  general,  is  not  con¬ 

vex)  .  In  the  p-version  of  the  finite  element  method  polynomials  are  used 

as  approximations  over  convex  subdomains.  This  critical  difference  gives 
the  p-version  a  more  rapid  rate  of  convergence  than  either  the  Ritz  method 
or  the  h-version. 

3.  The  h-p  version  of  the  finite  element  method.  In  this  approach  both  the 
degree  p  of  the  approximating  polynomial  and  the  maximum  diameter  h  of  the 
convex  subdomains  are  allowed  to  change. 
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The  p-version  of  the  finite  element  method  requires  families  of  polyno¬ 
mials  of  arbitrary  degree  p  defined  over  different  geometric  shapes.  Pol%’no- 
mials  defined  over  neighboring  elements  join  either  continuously  (are  in  C^) 
for  planar  or  three  dimensional  elasticity,  and  smoothly  (are  in  C^)  for  plate 
bending.  In  order  to  implement  the  p-version  efficiently  on  the  computer, 
these  families  should  have  the  property  that  computations  performed  for  an 
approximation  of  degree  p  are  re-usable  for  computations  performed  for  the 
next  approximation  of  degree  p  +  1.  We  call  families  possessing  this  property 
hierarchic  families  of  finite  elements. 

The  h-version  of  the  finite  element  method  has  been  the  subject  of  inten¬ 
sive  study  since  the  early  1950 's  and  perhaps  even  earlier.  Study  of  the 
p-version  of  the  finite  element  method,  on  the  other  hand,  began  at  Washington 
University  in  St.  Louis  in  the  early  1970's  and  led  to  a  more  recent  study  of 
tlie  h-p  version.  Researcli  in  the  p-version  (formerly  called  The  Constraint 
Method)  has  been  supported  in  part  of  the  Air  Force  Office  of  Scientific  Research 
since  1976. 


m 
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3 .  RESEARCH  ACCOI-IPLISHMENTS 

2.1.  Rate  of  convergence  of  the  p-version 

Extensive  computational  experiments  have  long  furnished  empirical  evidence 
that  the  rate  of  convergence  of  the  p-version  is  significantly  higher  than  that 
of  the  h-version  (See  for  example  [A,  5,  6]) 

The  following  two  theorems  provide  a  mathematical  explanation  for  the 
efficiency  of  the  p-version.  In  both  tlieorems  ..  is  a  bounded  polygonal  domain 
in  the  plane.  In  Theorem  1,  a  model  problem  for  the  case  is  considered,  and 
in  Theorem  2,  a  model  problem  for  the  case  is  considered.  In  both  cases, 
the  problems  are  singularity  problems,  that  is  the  smoothness  of  the  solutions 
are  governed  by  the  local  behavior  at  the  vertices  A.  of  the  polygons.  Suppose 
that  i.  is  a  angle  at  vertex  A.,  and  that  polar  coordinates  (r.,  :,)  are  used 

1  1  1  L 

at  A. .  The  solution  in  the  neighborhood  of  A.  is  of  the  form 
1  1 


Ci(r.)  e.c.) 


Y  . 


o^(r^)  »  g^dlog  r^l),  9  is  very 


smooth . 


Y^  depends  upon 


The  domain  Ti  is  assumed  to  be  triangulated  in  such  a  way  that  coincides 
with  a  vertex  of  a  triangle.  The  exact  assumptions  on  and  are  given 
in  [5]  and  [&]. 

Theorem  1  (model  problem  for  Case) 

-Au  +u=finQ  (1) 

U  =  0  01.  j.': 


Let  u  be  the  solution  to  (1)  in  the  weak,  sense,  and  let  u  be  the  finite 

P 

element  approximation  Co  u,  using  polynomials  of  degree  p  with  the  triangu¬ 
lation  S  fixed  (i.e.  u^  is  the  solution  to  (1)  using  the  p-version  of  the 
finite  element  method).  If  u  e  H  (H) ,  with  k  >  1,  then 


^  P 


where  c  >  0  is  arbitrary. 


u  =  min  (k-1,  2v^) ,  (2) 


Theorem  2  (model  problem  for  C  case) 


A“  w  =  f  in 


w  =  ■: —  =0  on  (clamped  edge) 


Let  w  be  the  solution  to  (2)  in  the  weak  sense,  and  let  w  be  the  solution  to 

P 

(3)  using  the  p-version  of  the  finite  element  method.  If  w  c  H  (2)  with 
k  >  2 ,  then 


!  |w  -  Wpl  I2  f  c  P  1  Iwj  Ij^,  u  “  min(k-2,  2(y^  -  D) 


Y^(a^)  •  smallest  positive  root  of  the  equation 


2  2  2 
sin  (y^  -  ^^“i  ~  ^^i  ~  =  0 


* 

Vi 
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The  error  estimates  in  (2)  and  (A)  can  be  compared  with  analogous  esti¬ 
mates  for  u,  and  w.  ,  the  solutions  by  the  h-version  of  the  finite  element 

h  n 

method.  Assume,  for  convenience  chat  in  (2)  k  -  1  ^  "^i’ 
k  -  2  2.  is  convergence  is  determined  by  the  nature  of  the 

singularities  at  comers.  Then  the  analogous  estimates  are 


!  iu  -  Uy,!  <  ch  Y  =  rain 

I  N  -  "hi  ^2, a  -  1  Ik’  ^  '^i 


(2') 

(A’) 


1/2  —1/2 

If  N  is  the  number  of  degrees  of  freedom  then  p-N  “.h-N  and 

11“  -  “pi  I i.n-  <>(«■’'>.  11“  -  “hlli.s- 


1“  -  “p!  l2,n  '  l!“  -  “Jl2,a  - 


Therefore,  if  the  criterion  used  to  compare  methods  is  the  number  of  degrees 
of  freedom  required  to  achieve  a  given  error  in  energy,  then  the  rate  of  con¬ 
vergence  of  the  p-version  is  twice  that  of  the  h-version.  This  result  provides 
a  rigorous  proof  for  the  extensive  computational  evidence  that  has  been  gathered, 
and  explains  in  part  the  efficiency  of  the  p-version. 

Also,  it  has  been  shown  in  [7,  8]  that  for  the  h-p  version  under  the  assump¬ 
tions  in  Theorems  1  and  2  if  optimal  mesh  refinement  (not  necessarily  quasi 
uniform)  is  combined  with  optimal  p-distribution  then  the  estimate  becomes  (in 
the  C  case) 


wnere  t,-^  are  constants  and  ••  _  1/3.  Therefore,  for  the  combined  h-p  version 
the  asymptotic  rate  of  convergence  becomes  exponential.  It  is  demonstrated  in 
[7]  that  an  exponential  rate  of  convergence  is  achieved  for  the  p-version  in 
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c:':e  case  cr  problems  with  smooth  solutions  and  that  problems  with  stress  singu¬ 
larities  can  be  made  to  behave  almost  as  problems  with  smooth  solutions  providL-a 
that  properly  graded  meshes  are  used  for  the  p-version. 

Figure  1  illustrates  graphically  the  rates  of  convergence  of  the  h-,  p- 
and  h-p  versions. 

2.2.  Hierarchic  Families  of  Solid  Finite  Elements 

In  order  to  implement  the  p-version  efficiently,  families  of  finite  elements 
are  needed  with  the  hierarchic  property:  computations  performed  for  an  approxi¬ 
mation  of  order  p  should  be  re-usable  when  raising  the  order  to  p  +  1.  More 
specif icially ,  the  stiffness  matrix  corresponding  to  the  polynomial  approxi¬ 
mation  of  degree  p  should  be  a  submatrix  of  the  polynomial  approximation  of 
degree  p  +  1.  In  terms  of  basis  functions,  this  implies  that  the  basis  func¬ 
tions  for  a  pth  order  approximation  should  be  a  subset  of  the  basis  functions 
for  a  (p  +  l)st  order  approximation. 

Hierarchic  families  for  triangles  both  in  the  case  and  in  the  case 
are  described  in  detail  in  [9,  10,  11,  12,  13]  We  now  briefly  describe  some 
hierarchic  families  of  polynomials  for  various  three  dimensional  shapes.  All  of 
these  families  are  globally  in  C*^. 

tetrahedron.  A  hierarchic  family  for  the  tetrahedron  can  be  constructed  from 
the  hierarchic  family  for  the  triangle  by  using  natural  coordinates.  The 
linear  element  is: 


triangular  prism  A  hierarchic  family  for  the  triangular  prism  can  be 
constructed  by  combining  the  hierarchic  families  for  the  triangle  (in  natural 
coordinates)  and  for  the  square  (in  rectangular  coordinates). 

Details  are  given  in  [l4]  and  [is]- 


square  pyramid.  It  can  be  shown  [15]  that  no  hierarchic  family  exists 
which  consists  of  polynomials  alone.  However  it  is  possible  to  supplement 
certain  rational  functions  in  such  a  way  that  the  element  of  degree  p  contains 
a  complete  polynomial  of  degree  p  and  additional  rational 
functions.  Furthermore,  because  of  the  special  form  of 
these  rational  functions,  integration  of  all  shape  func¬ 
tions  which  appear  in  the  elemental  stiffness  matrix  can 
be  performed  in  closed  form.  No  numerical  quadrature 
is  required.  See  [l5]  for  further  details. 

Combinations 

Hierarchic  families* of  the  different  shapes  shown  here  have  been  con¬ 


(&ai) 


structed  so  that  they  join  together  continuously.  Thus  geometries  such  as  the 
ones  shown  below  can  be  modeled  using  as  few  elements  as  possible.  In  this 
way  the  subdivision  of  the  polyhedral  domain  T.  can  be  made  verv  coarse  and 


accuracy  can  be  obtained  by  increasing  p,'  without  using  added  degrees  of 
freedom  to  describe  the  geometry.  As  has  been  shown  (at  least  for  two  diinen 
sional  problems)  the  rate  of  convergence  in  the  p-version  is  twice  that 
obtained  for  the  h-version.  Therefore  subdividing 


in  this  way  leads  to  computational  efficiency,  since  it  makes  maximum  use  of 
the  p-version.  It  also  requires  significantly  simpler  input. 


2.3  t)uanticies  of  Special  Interest  -  Extraction  Techniques 


Often  the  main  purpose  of  a  finite  element  analysis  is  to  obtain  values  of 
a  few  important  quantities  with  a  high  accuracy.  In  structural  mechanics, 
for  example,  the  values  of  displacements  or  stresses  in  a  small  number  of 
designated  areas,  or  the  stress  intensitv  factor  at  a  small  number  of  points 
is  of  critical  importance  for  design.  Both  the  h-  and  p-  versions  of  the 
finite  element  method  give  approximations  for  these  values.  However,  it  is 
much  more  efficient  to  use  a  post-processing  technique  which  uses  weighted 
averages  of  values  taken  directlv  from  the  finite  element  approximation.  The 
post-processor  determines  these  quantities  of  special  interest  much  more 
accurately.  In  particular,  ii^hen  stresses  are  computed  pointwise  as  derivatives 
of  displacements  in  the  p-version,  they  may  exhibit  (some  times  severe)  oscil¬ 
latory  behavior.  In  the  case  of  the  centrally  cracked  panel  shown  in  Figure 
2(h) ,  for  example  the  normal  stresses  along  the  x-axis  computed  for  poly¬ 

nomial  orders  ranging  from  1  to  7  are  shown.  The  oscillatory  behavior  near 


the  crack  tip  singularity  is  evident.  Two  techniques  are  being  developed  tor 


tile  nos tp rocess ing  of  quantities  of  special  interest.  Both  teciiniques  are 
ve rp'  accurate,  vielding  approxiniations  that  are  of  the  same  order  of  accurac’. 
as  the  strain  energy.  (This  is  the  square  of  the  error  in  energy  norm.) 

2.3.1  Use  of  Functional  Forms 

This  approach  is  based  on  the  idea  that  the  functional  forms  of  the  quan 
titles  of  interest  and  generally  known.  In  the  case  of  a  centrallv  cracked 
panel,  for  example,  the  displacements  u(z),  v(z),  and  the  stresses  '  (z), 

T  (z)ando(z)  (z=x+iv)  are  given  in  terms  of  two  functions. 

XV  V 


;  (z)  =  ;  (z)  + 


4'^(z) 


/r 


b  .  z-^  + 


J-  1/2 


(z)  =- ;  (z)  + 


^  b.2-*  +  -  a.zJ  - 

°  J7  >  l-O  1 


who  re 


i  (z)  =  r:  b.z“ 

°  j=0  J 


, (z)  =  2  a.z“ 

‘  j-O  J 


are  holomorphic  functions.  Approximate  values  for  the  coefficients  a  ,  b 

1  j 

are  determined  using  the  displacements  u  (z),  v  (z)  computed  by  the  p-version 

P  P 

of  the  finite  element  method. 

This  technique  was  used  in  [l6]  to  obtain  improved  values  of  o  ^  along 
tlie  x-axis  in  the  case  of  a  centrally  cracked  panel. 


A  general  form  for  post  processing  calculations  is  given  in  [  17,  18  ]. 
Green's  theorem  and  generalized  influence  functions  are  used  together  with 
smooth  cut-off  functions  and  blending  functions  in  order  to  calculate  higher 
derivatives  of  the  unknown  function  and  also  stress  intensity  factors. 

These  techniques  have  been  applied  to  post  processing  of  the  u  and  v 

P  P 

displacement  fields  obtained  from  the  p-version.  Figure  3  shows  two  meshes 
for  the  p-version  finite  element  analysis  of  a  cracked  plate.  Figures  4  and 
5  compare  the  absolute  error  in  computed  directly  from  the  p-version  data 
with  that  obtained  by  the  extraction  method.  The  stress  itself  is  infinite 
at  the  crack  tip  (although  the  strain  energy  is  finite).  The  extraction 
method  shown  in  B,  however,  surpresses  the  spurious  behavior  of  away  from 
the  crack  tip. 


2.4  A  Sample  Problem:  Analysis  of  an  Attachment  Lug 


An  example  of  the  use  of  tlie  p-version  for  estimating  and  controlling 


the  error  in  the  approximation  is  the  analysis  of  an  attachment  lug,  shown 
below.  This  problem  was  eonsidcced  in  detail  in  [  19 ,  20]  and  is  of  considerable 


interest  in  the  aerospace  industry  since  similar  structural  details  are  used 


Attachment  lug 


The  dimensions  shown  are  in  inches.  The  thickness  of  the  lug  is  0.5 
inches  (12.7mmj.  The  material  is  assumed  to  be  isotropic,  elastic,  with  mod¬ 
ulus  of  elasticity  of  30,000  ksi  (207000  MPa)  and  Poisson's  ratio  of  0.3. 

Plane  stress  conditions  are  assumed.  The  lug  is  fixed  along  A-B  and  is  loaded 
by  a  tightly  fitting  circular  pin  through  the  hole.  The  pin  imposes  a  load 
of  10.0  kips  (44.5  kN).  The  line  of  action  of  the  applied  load  passes  through 
the  center  of  the  circular  hole  and  is  inclined  to  the  horizontal  by  45  decrees 
This  force  is  to  be  represented  by  a  sinusoidal  normal  pressure  distribution 
acting  on  the  inside  surface  of  the  hole.  The  pressure  distribution  is  to 
act  within  the  range  of  +90.0  degrees  to  -90.0  degrees  from  the  direction  of 
the  applied  load,  as  shown  below.  The  line  force  q  is  therefore: 

q  =  I’sin  0  <  <  •>t 

Where  Q  =  5.093  kips/in  (891  kN/m)  and  •  is  measured  from  the  axis  counter¬ 
clockwise  as  shown  below.  The  goal  of  the  computation  is  to  determine  the 
ma:\imum  principal  stress  in  the  neighborhood  of  the  hole. 


Sinusoidal  loading 


In  the  p-version  of  the  finite  element  method  the  mesh  is  normally  chosen 
to  represent  the  plane  elastic  bodv  with  a  few  elements  as  possible.  An  impor¬ 
tant  exception  is  the  neighborhood  of  singular  points  which  should  be  isolated 
by  one  or  more  layers  of  small  elements.  In  this  problem  points  A  and  B  are 
singular  points.  Under  the  assumptions  of  linear  elasticity,  the  stress  is 
infinity  at  these  points.  According  to  our  problem  statement,  points  A  and  B 
are  outside  of  the  region  of  the  interest,  however.  Otherwise  it  would  have 
been  necessary  to  define  either  the  gcjal  of  the  computation  or  the  support 
conditions  differently.  The  purpose  of  the  small  elements  A  and  B  is  to  reduce 
the  spreading  of  approximation  error  caused  by  the  singular  points.  The  spread¬ 
ing  of  error  caused  by  singular  points  is  called  'pollution'.  Our  mesh  design 
is  shown  in  Fig.  6.  The  mesh  used  in  [20 3  when  the  h-version  is  employed  is 
shown  in  Figure  7.  In  [l9j  techniques  wliich  assess  the  quality  of  finite 
element  solutions  are  performed  through  the  follo\;ing  tests; 

(a)  Estimate  the  relative  error  in  energy  norm; 

(b)  Perform  an  overall  equilibrium  test; 

(c)  Perform  an  element  by  element  equilibrium  test; 

(d)  Perform  a  convergence  test  on  the  magnitude  and  orientation 
of  principal  stress  components  of  node  18; 

(e)  Perform  a  convergence  test  on  the  location,  magnitude  and 
orientation  of  the  maximum  principal  stress. 

Here  we  give  the  results  for  only  (a)  and  (e) . 

(a)  Estimated  relative  error  in  energy  norm. 

The  number  of  degrees  of  freedom,  N,  the  computed  strain  energy,  U(u  ),  the 

— FL 

estimated  exact  strain  energy,  the  rate  of  convergence  S  and  the  estimated 
relative  error  are  shown  in  Table  1. 


Table  1 


t  inate  i'll’  error  in  ener;^;'  norm  based  on  p-extens  ion 


p_ 

N 

Computed  Strain 
Energy* 

Estimated  Exact 
Strain  Energy* 

Percent 
Error  E: 

1 

36 

0.173592D-01 

- 

- 

- 

2 

100 

0.276935D-01 

- 

- 

- 

3 

170 

0. 283543D-01 

0.285992D-01 

2.47 

9.25 

4 

266 

0.285719D-01 

0. 28734 3D-01 

1.90 

7.52 

5 

388 

0. 286491D-01 

0.287139D-01 

2.08 

4.75 

ft 

5  36 

0. 286839D-01 

0.287277D-O1 

1.81 

3.91 

7 

710 

0.286906D-01 

0.2869290-01 

4.90 

0.89 

8 

910 

0.286928D-01 

0. 2869420-01 

3.  74 

0.71 

*  inch-kip  units 

It  is  seen  that  the  relative  error  in  energy  norm  is  '  timated  to  be  under 
one  percent  for  p  =  7  and  8. 

(e)  The  maximum  principal  stress 

In  elastic  stress  analvsis  the  "oal  of  finite  elemi^nt  I'omputation  is  usually 
to  find  the  location  and  magnitude  of  the  maximum  principal  stress,  or  the 
maximum  shear  stress.  In  the  p-version  this  is  done  by  computing  tlie  stress 
components  on  a  fine  grid  and  sorting  them  either  according  to  the-  maximum 
principal  stress  or  the  maximum  shear  stress.  The  grids  are  imposed  on  the 
standard  elements  and  mapped  onto  the  'real'  element  by  the  mapping  function 
of  the  element. 

In  the  case  of  the  lug  problem  we  chose  a  6  x  6  grid  per  element  and 
searched  for  the  maximum  principal  stress  on  elements  3  to  10.  The 
results  are  shown  in  Table  2.  Convergence  is  clearly  visible.  To  find  a 
sharper  definition  of  the  maximum  principal  stress  we  also  specified  a  10  x  If 
grid  on  element  4.  The  maximum  principal  stress  on  this  grid  was  found  to  be 
24.91  ksi,  its  orientation  41.1  degrees  from  the  x^  axis  and  its  location: 
x^  =  -0.812  in;  =  0.951  in. 
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Table  2 

Location  and  Magnitude  of  the  Maximum 
i’rincipal  Stress  (6  x  6  ”rid  per  element) 


£ 

Element 

No. 

Location 

Principal  Stresses  (ksi) 

CT  T 

1  2 

1 

4 

-0.884 

0.  884 

21.67 

3.  71 

43.9 

2 

4 

-0.625 

1.083 

20. 88 

1.99 

35.0 

3 

4 

-0.625 

1.083 

23.  73 

1.84 

35.0 

4 

4 

-0.625 

1.083 

23.64 

0.58 

35.2 

5 

<4 

-0. 761 

0.Md2 

23.52 

-0.58 

39.5 

6 

4 

-0.761 

0.992 

24.46 

-0.77 

39. 2 

7 

4 

-0.761 

0.992 

24.78 

-0.80 

38.  7 

8 

4 

-0. 761 

0.992 

24.87 

-0.04 

38.4 

Orientat 

ion  of  first 

principal 

axis  relative 

to  axis 

x^ ,  counterc 

lockwise 

in  de;4rees. 

2 . 5  Computer  Implementations 

Our  earlier  work  on  researcli  and  develonment  of  the  n-version  have  resulted 
in  several  computer  implementations.  As  our  research  progresses,  new  procedures 
and  ideas  are  continually  incorporated  into  these  codes  making  them  more  effic¬ 
ient.  Two  of  the  codes  (FIESTA/3D  and  PROBE)  are  in  commercial  use. 

2.5.1  COMET-X 

COMET-X  is  an  experimental  computer  code  which  implements  the  p-version 
of  the  finite  element  method  by  using  the  hierarchic  families  which  have  been 
constructed.  COMET-X  is  maintained  by  the  Center  for  Computational  Mechanics 
at  Washington  University.  COMET-X  can  be  used  as  a  code  to  implement  the  h- 
version  as  well  simply  by  fixing  the  polynomial  order  p  and  refining  the  mesh. 

COMET-X  currently  has  the  following  capabilities; 

A.  Element  types:  Stiffeners,  triangular  elements,  triangular  elements  with 

one  side  curved,  rectangular  elements,  solid  elements  of  the  shapes  described 
earlier. 


B.  Types  of  Analvsis:  Laplace  and  Poisson  equations,  plane  elasticity, 
tenpetature  distribution  in  l-dimens ions . 

C.  Special  Capabilities:  non-uniform  p-dis tribution ,  elastic  fracture  mech¬ 
anics  computations  in  two  dimensions,  nearly  incompressible  solids,  linear 
boundary  layer  problems. 

D.  Pre-  and  Post  processing  capabilities  including  graphics  and  visual  display 
Tlie  capabilities  and  usage  of  COMET-X  are  described  in  detail  in  [2l]. 

2.5.2  FIESTA/3D 

FIESTA/3D  is  a  software  system  for  static  analysis  of  solid  structures 
based  on  the  p-version  of  the  finite  element  method.  It  is  marketed  by 
McDonnell-Douglas  Automation  Company.  Some  of  the  advanced  features  available 
on  FIESTA/ 3D  are: 

•  controllability  of  the  quality  of  the  solution 

•  efficient  modeling 

•  treatment  of  stress  singularities 

•  a  method  for  surface  identification 

2.5.3  PROBE 

PROBE  is  an  advanced  computer  implementation  of  the  p-version  for  2-dimen- 
sional  analysis  currently  under  development  at  NOETIC  Technologies  in  St.  Louis 
Missouri.  Some  of  PROBE's  unique  features  are: 

•  automatic  error  estimators.  These  estimators  provide  initial  feedback 
on  solution  quality  by  computing  element-by-element  and  edge-by-edge 
equilibrium  checks. 

•  interactive  error  estimators.  These  estimates  offer  immediate  feed¬ 
back  on  solution  quality  by  computing  action/ reaction  checks  overall 
equilibrium  checks  and  a  convergence  trajectory  of  stresses  and  strains 


•  load  tracking,  for  extracting  free-body  diagrams 


•  pinpoint  solutions  which  give  specific  results  anywhere  in  the  model 

•  Mode  I  (Symmetric)  and  Mode  II  (antisymmetric)  Stress  Intensity  Factors 
for  Fracture  Mechanics 

•  Precise  curve  definitions  which  eliminate  mesh  refinement  near  cutouts 
and  provide  more  precise  solutions  in  critical  regions. 


(The  relative  error  values  shown  are  typical  for  certain  engineering  proble 
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Zleriena  Method",  Svr'.oosi'art  or.  Ao'/a.r'.oe-  ar.c  r'e:'...5  in  r 1 1 ..  1 

ana  Solid  Mecnkr.ics,  October  0-7,  1982,  '■.'as'nin^ ton ,  1.  w. 

"?-ConvergenC  Polynomial  Approximations  in  by  Douulas  hang  an 

I.  Norman  K.aC2,  Fourth  Texas  Symposium  on  Approximation  Theory.  Denartme 
of  .Mathematics,  Texas  A  i  M  Universtiy,  College  Station  Texas  77S-.3, 
January  17-21,  1983. 

"Design  Aspects  of  Adaptive  Finite  Element  Codes"  by  D.  W.  Vanq,  I.  N. 
Katz  and  M.  Z.  Qian,  ASCE-EMD  (American  Society  of  Civil  Eneineers- 
Engineering  Mechanics  Division)  Soecialitv  Conference,  '^u’^due  ve’-  =  - -- 

May  25-28,  1983.  ‘  ‘  ’  ’ 

"Smoothing  Stresses  Computed  Pointwise  by  the  p-version  of  the  Finite 
Element  Method"  by  I.  Norman  Katz  and  Xing-ren  Ying,  SIAM  1983  National 
Meeting,  Denver,  Colorado,  June  6-8,  1983. 

"T!ie  Use  or  I'rder  Po  L  vni'mi  a  Is  in  the  Numerical  Solution  o;  Pirti.. 

Differencial  Kquations",  a  Mini  Symposium.  1.  Norman  KaCr..  a  rc  in  i r  o--. 
Chairman;  "The  h-p  version  of  the  Finite  Element  Method",  I.  Babuska, 

B,  Szabo,  K.  Izadpanah,  W.  Gui,  and  B.  Guo;  "A  Pseudospectral  Legendre 
Method  for  Hyperbolic  Equations",  D.  Gottlieb  and  H.  Tal-Ezer;  "The 
Approximation  Theory  for  the  p-version  of  the  Finite  Element  .Method", 

Milo  Dorr;  "On  the  Robustness  of  Higher  Order  Elements",  M.  Vogelius; 

SI.-VM  Summer  Meeting,  L'niver.sitv  of  Washington,  Seattle,  Washington, 

July  16-20,  1984. 

"Implementation  of  a  Finite  Element  Software  System  with  h-  and  p-exten- 
sion  capabilities'  ,  by  Barna  A.  Szabo,  to  be  presented  at  the  8th  Invita¬ 
tional  Symposium  on  the  Unification  of  Finite  Element-Finite  Differences 
and  the  Calculus  of  Variations,  University  of  Connecticut,  Storrs,  Connecti 
cut,  .May  3,  1985. 

"Stress  Singularities  at  Angular  Corners  of  Composite  Plates",  Xing-Ren  Yin 
and  1.  Norman  Katz,  SIAM  Spring  Meeting,  Pittsburgh,  PA,  June  1985. 

"An  Overview  of  the  p-version  of  the  Finite  Element  Method,"  by  I,  Norman 
Katz,  invited  colloquium  presentation  at  the  Department  of  Mathematics, 

West  Virginia  University,  Morgantown,  West  Virginia,  June  27,  1985. 


I  5.3  Visits  and  Seminars  at  Government  Laboratories 

I  1.  "Advanced  Stress  Analysis  Technology"  3Y  3.  A.  Szabo  and  I.  X.  Katz, 

presented  on  September  8,  1977  at  the  'Air  Force  Flight  Dynamic  Labora¬ 
tory,  Wright-Patterson  Air  Force  Base. 

abstract 

With  one  exception,  all  finite  element  software  systems  have  element 
j  libraries  in  which  the  approximation  properties  of  elements  are  frozen. 

,  The  user  controls  only  the  number  and  distribution  of  finite  elements. 

I  The  e.xception  is  an  experimental  software  system.,  developed  at  '-.'ashington 

I  Tniversity .  This  system,  called  COMET-M,  employs  conforming  elements 

based  on  complete  polynomials  of  arbitrary  order.  The  elements  are  hierar¬ 
chic,  i.e.  the  stiffness  matrix  of  each  element  is  em.bedded  in  the  stiff¬ 
ness  matrices  of  all  higher  order  elements  of  the  same  kind.  The  user  con¬ 
trols  not  only  the  number  and  distribution  of  finite  elements  but  their 
approximation  properties  as  well.  Thus  convergence  can  be  achieved  on 
fixed  mesh.  This  provides  for  very  efficient  and  highly  accurate  approxi¬ 
mation  and  a  new  method  for  computing  stress  intensity  factors  in  linear 
elastic  fracture  mechanics.  The  theoretical  developments  are  outlined, 
numerical  examples  are  given  and  the  concept  of  an  advanced  self -adant ive 
finite  element  software  system  is  presented. 

2.  "The  Constraint  Method  for  Finite  Element  Stress  .Analysis",  bv  I.  N. 

Katz,  presented  at  the  National  Bureau  of  Standards,  Apnlied  Mathe¬ 
matics  Divison  on  October  19,  1977. 

abstract 

In  conventional  approaches  to  finite  element  stress  analysis  accuracy 
is  obtained  by  fixing  the  degree  p  of  the  approximating  polynom.ial  and  by 
allowing  the  maximum  diameter  h  of  elements  in  the  triangulation  to  approach 
zero.  An  alternate  approach  is  to  fix  the  triangulation  and  to  increase  the 
degrees  of  approximating  polynomials  in  those  elements  where  more  accuracy 
is  required.  In  order  to  implement  the  second  approach  efficiently  it  is 
necessary  to  have  a  family  of  finite  elem.ents  of  arbitrary  polynomial  degree 
p  with  the  property  that  as  much  information  as  possible  can  be  retained  from 
the  pth  degree  approximation  when  computing  the  (p+l)st  degree  approximation. 
Such  a  HIEPwARCKIC  family  has  been  formulated  with  p  ^  2  for  problems  in  plane 
stress  analysis  and  with  p  ^  5  for  problems  in  plate  bending.  The  family  is 
described  and  numerical  examples  are  presented  which  illustrate  the  efficiency 
of  the  new  method. 

3.  "The  p-version  of  the  Finite  Element  Method",  by  I.  N.  Katz  and  3.  A. 
Szabo,  to  be  presented  at  Air-Force  Flight  Dynamics  Laboratory  Wright- 
Patterson  Air  Force  Base  on  April  23,  1981  (tentative  date) . 

abstract 


The  theoretical  basis  of  the  p-version  of  the  finite  element  method  has 
been  established  only  quite  recently.  Nevertheless,  the  p-version  is  already 


seer,  zo  be  the  most  promising  approach  for  implementing  adaptivity  in  practi¬ 
cal  computations.  The  main  theorems  establishing  asymptotic  rates  of  conver¬ 
gence  for  the  p-version,  some  aspects  of  the  algorithmic  structure  of  p-version 
computer  codes,  numerical  experience  and  a  posteriori  error  estimation  will  be 
discussed  from  the  mathematical  and  engineering  points  of  view. 

.Attended  the  meeting  sponsored  by  AFOSR  on  "The  Impact  of  Large  Scale 
Computing  on  Air  Force  Research  and  Development",  at  Kirtland  Air  Force 
Base,  Albuquerque,  New  Mexico,  April  4-6,  1984. 


